Install/load necessary R packages
Set date range
min_year = 2017
max_year = 2021
More information about peatland and upland water elevation data collected at the Marcell Experimental Forest can be found on the Environmental Data Initiative Repository Portal.
Read in the most up-to-date data from the following Box file
path:
External-MEF_DATA/Hydro/WaterTableElevation/Bogwells/L1
External-MEF_DATA/EDI/peatland_water_table_daily/edi.562.2/data_objects/Clean the data into a tidy format
#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/Hydro/WaterTableElevation/Bogwells/L1"
#read in the raw data
bogwell <- read_csv(here(filepath, "L1_Peatland_well_daily_history_2021forMia.csv"))
#clean the data
##note: we do not want any E values for the flag column
bogwell_clean <- bogwell %>%
clean_names() %>%
mutate(year = format(as.POSIXct(date, format = '%Y-%m-%d',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year)) %>%
subset(year >= min_year & year <= max_year)
#save clean data CSV to use in future RMD files
write.csv(x = bogwell_clean,
file = file.path(here("intermediate_data", "bogwell_clean.csv")),
row.names = FALSE)
Plot all sites together
#plot
p_bogwell <- ggplot(data = bogwell_clean) +
geom_line(aes(x = date, y = wte, color = peatland)) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")"),
color = "Site"
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#static plot
p_bogwell
bogwell_clean_s1 <- bogwell_clean %>%
subset(peatland == "S1")
#plot
p_bogwell_s1 <- ggplot(data = bogwell_clean_s1) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S1 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s1,
tooltip = "text")
bogwell_clean_s2 <- bogwell_clean %>%
subset(peatland == "S2")
#plot
p_bogwell_s2 <- ggplot(data = bogwell_clean_s2) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S2 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s2,
tooltip = "text")
bogwell_clean_s3 <- bogwell_clean %>%
subset(peatland == "S3")
#plot
p_bogwell_s3 <- ggplot(data = bogwell_clean_s3) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S3 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s3,
tooltip = "text")
bogwell_clean_s4 <- bogwell_clean %>%
subset(peatland == "S4")
#plot
p_bogwell_s4 <- ggplot(data = bogwell_clean_s4) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S4 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s4,
tooltip = "text")
bogwell_clean_s5 <- bogwell_clean %>%
subset(peatland == "S5")
#plot
p_bogwell_s5 <- ggplot(data = bogwell_clean_s5) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S5 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s5,
tooltip = "text")
bogwell_clean_s6 <- bogwell_clean %>%
subset(peatland == "S6")
#plot
p_bogwell_s6 <- ggplot(data = bogwell_clean_s6) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S6 Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_s6,
tooltip = "text")
bogwell_clean_boglk <- bogwell_clean %>%
subset(peatland == "BOGLK")
#plot
p_bogwell_boglk <- ggplot(data = bogwell_clean_boglk) +
geom_line(aes(x = date,
y = wte,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", wte))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("BOGLK Daily Peatland Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_bogwell_boglk,
tooltip = "text")
Function for loading Campbell DataLogger .dat files:
readCDL = function(file){
# read data file starting on 5th line
dat <- read.csv(file, sep=",",header=FALSE,skip=4,na.strings="NAN",stringsAsFactors=F)
# Read in just the header line (l2)
# unlist the line, and remove quotes
h <- readLines(file, n=2)[2]
n <- as.factor(unlist(strsplit(h, ",")) )
n2 <- gsub('"', "", n)
# assign column names to dataframe
colnames(dat) = n2
return(dat)
}
Constants to convert between feet and meters:
MetersPerFoot = 0.3048
FeetPerMeter = 3.28048
Read in the 2022 data from the following Box file path:
External-MEF_DATA/DataDump/RealTimeData
Clean the 2022 data into a tidy format
#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/RealTimeData"
#read in the raw data
##note: this data is only 2022
s1_se <- readCDL(file = here(filepath, "S1-EM3_Table1.dat"))
#clean the data
s1_se_clean <- s1_se %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
wt_elev_ft = as.numeric(wt_elev_ft),
wt_elev_m = wt_elev_ft * MetersPerFoot)
External-MEF_DATA/DataDump/annual_appende_logger_files/2019#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
#set the year
year = "2019"
#set the file name
file = "S1-EM3_Table1.dat"
#read in 2019 data
s1_2019 <- readCDL(file = here(filepath, year, file))
#clean 2019 data
s1_2019_clean <- s1_2019 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
wt_elev_ft = as.numeric(wt_elev_ft),
wt_elev_m = wt_elev_ft * MetersPerFoot)
External-MEF_DATA/DataDump/annual_appende_logger_files/2020#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
#set the year
year = "2020"
#set the file name
file = "S1-EM3_Table1.dat"
#read in 2020 data
#note: this data is 2020 - 2022
s1_2020 <- readCDL(file = here(filepath, year, file))
#clean 2020 data
s1_2020_clean <- s1_2020 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
wt_elev_ft = as.numeric(wt_elev_ft),
wt_elev_m = wt_elev_ft * MetersPerFoot)
RealTimeData folder
and annual_appnded_logger_files folder in Boxs1_all <- rbind(s1_se_clean, s1_2019_clean, s1_2020_clean) %>%
#remove NA values
filter(!is.na(wt_elev_m)) %>%
#extract year
mutate(year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year)) %>%
#remove duplicate entries
unique()
#set min and max year
min_year = min(s1_all$year)
max_year = max(s1_all$year)
#plot
p_s1 <- ggplot(data = s1_all) +
geom_line(aes(x = datetime,
y = wt_elev_m,
group = 1,
text = paste0("Date: ", datetime, "\n",
"Water Table Elevation (m): ", round(wt_elev_m, digits = 3)))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S1 Shaft Encoder Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_s1,
tooltip = "text")
External-MEF_DATA/DataDump/RealTimeData
MaxWTElev column#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/RealTimeData"
#read in the raw data
s2_se <- readCDL(file = here(filepath, "S2_bogwell_S2BW_Daily.dat"))
#clean the data
s2_se_clean <- s2_se %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = max_wt_elev_ft * MetersPerFoot,
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
#set min and max year
min_year = min(s2_se_clean$year)
max_year = max(s2_se_clean$year)
#plot
p_s2 <- ggplot(data = s2_se_clean) +
geom_line(aes(x = datetime,
y = max_wt_elev_m,
group = 1,
text = paste0("Date: ", datetime, "\n",
"Water Table Elevation (m): ", max_wt_elev_m))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m)",
title = paste0("S2 Bogwell Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
#interactive plot
ggplotly(p_s2,
tooltip = "text")
S3 was omitted due to unreliable data.
External-MEF_DATA/DataDump/annual_appende_logger_files/2020#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2020"
file = "S4_bogwell_S4BW_Daily.dat"
#read in 2020 data
s4_2020 <- readCDL(file = here(filepath, year, file))
#clean the data
s4_2020_clean <- s4_2020 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
#set min and max year
min_year = min(s4_2020_clean$year)
max_year = max(s4_2020_clean$year)
#plot
p_s4 <- ggplot(data = s4_2020_clean) +
geom_line(aes(x = datetime,
y = max_wt_elev_m,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", max_wt_elev_m))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m)",
title = paste0("S4 Bogwell Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
#interactive plot
ggplotly(p_s4,
tooltip = "text")
External-MEF_DATA/DataDump/annual_appende_logger_files/2020#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2020"
file = "S5_bogwell_S5BW_Daily.dat"
#read in 2020-2022 data
s5_2020 <- readCDL(file = here(filepath, year, file))
#clean the data
s5_2020_clean <- s5_2020 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
Read in the 2022 data from the following Box file path:
External-MEF_DATA/DataDump/RealTimeData
Clean the 2022 data into a tidy format
#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/RealTimeData"
#read in the raw data
##note: this data is 2021 - 2022
s5_se <- readCDL(file = here(filepath, "S5_bogwell_S5BW_Daily.dat"))
#clean the data
s5_se_clean <- s5_se %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
RealTimeData folder
and annual_appnded_logger_files folder in Boxs5_all <- rbind(s5_2020_clean, s5_se_clean) %>%
#remove NA values
filter(!is.na(max_wt_elev_m)) %>%
#remove duplicate observations
distinct()
#set min and max year
min_year = min(s5_se_clean$year)
max_year = max(s5_se_clean$year)
#plot
p_s5 <- ggplot(data = s5_all) +
geom_line(aes(x = datetime,
y = max_wt_elev_m,
group = 1,
text = paste0("Date: ", date, "\n",
"Water Table Elevation (m): ", max_wt_elev_m))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m)",
title = paste0("S5 Bogwell Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
#interactive plot
ggplotly(p_s5,
tooltip = "text")
Read in the 2022 data from the following Box file path:
External-MEF_DATA/DataDump/RealTimeData
Clean the 2022 data into a tidy format
#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/RealTimeData"
#read in the raw data
##note: this data is 2021 - 2022
s6_se <- readCDL(file = here(filepath, "S6_bogwell_S6BW_Daily.dat"))
#clean the data
s6_se_clean <- s6_se %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
External-MEF_DATA/DataDump/annual_appende_logger_files/2019#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2019"
file = "S6_bogwell_S6BW_Daily.dat"
#read in 2019 data
s6_2019 <- readCDL(file = here(filepath, year, file))
#clean 2019 data
s6_2019_clean <- s6_2019 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m))
External-MEF_DATA/DataDump/annual_appende_logger_files/2020#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2020"
file = "S6_bogwell_S6BW_Daily.dat"
#read in 2020-2022 data
s6_2020 <- readCDL(file = here(filepath, year, file))
#clean 2020-2022 data
s6_2020_clean <- s6_2020 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
max_wt_elev,
max_wt_elev_m) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = max_wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = as.numeric(max_wt_elev_m))
RealTimeData folder
and annual_appnded_logger_files folder in Boxs6_all <- rbind(s6_se_clean, s6_2019_clean, s6_2020_clean) %>%
#remove one row where WTE = 0
subset(max_wt_elev_m > 0) %>%
#remove NA values
filter(!is.na(max_wt_elev_m))
#set min and max year
min_year = min(s6_all$year)
max_year = max(s6_all$year)
#plot
p_s6 <- ggplot(data = s6_all) +
geom_line(aes(x = datetime,
y = max_wt_elev_m,
group = 1,
text = paste0("Date: ", date, "\n",
"Max Water Table Elevation (m): ", round(max_wt_elev_m, digits = 3)))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("S6 Shaft Encoder Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_s6,
tooltip = "text")
Read in the 2020-2022 data from the following Box file path:
External-MEF_DATA/DataDump/RealTimeData
Clean the 2020-2022 data into a tidy format
#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/RealTimeData"
#read in the raw shaft encoder data
boglk_se <- readCDL(file = here(filepath, "BLF_met_BogLakeW.dat"))
#clean the data
boglk_se_clean <- boglk_se %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = max_wt_elev_ft * MetersPerFoot,
max_wt_elev_m = as.numeric(max_wt_elev_m),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year))
External-MEF_DATA/DataDump/annual_appende_logger_files/2019#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2019"
file = "BLF_met_BogLakeW.dat"
#read in 2019 data
boglk_2019 <- readCDL(file = here(filepath, year, file))
#clean 2019 data
boglk_2019_clean <- boglk_2019 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = max_wt_elev_ft * MetersPerFoot,
max_wt_elev_m = as.numeric(max_wt_elev_m))
External-MEF_DATA/DataDump/annual_appende_logger_files/2020#set file path to read in data from Box
filepath = "/Users/miaforsline/Library/CloudStorage/Box-Box/External-MEF_DATA/DataDump/annual_appended_logger_files"
year = "2020"
file = "BLF_met_BogLakeW.dat"
#read in 2020 data
boglk_2020 <- readCDL(file = here(filepath, year, file))
#clean 2020 data
boglk_2020_clean <- boglk_2020 %>%
clean_names() %>%
#select relevant columns
select(timestamp,
wt_elev) %>%
#rename columns
rename(datetime = timestamp,
max_wt_elev_ft = wt_elev) %>%
#reformat columns into appropriate classes
mutate(datetime = as.POSIXct(datetime, format = '%Y-%m-%d %H:%M:%S',
tz = "GMT"),
date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
max_wt_elev_ft = as.numeric(max_wt_elev_ft),
max_wt_elev_m = max_wt_elev_ft * MetersPerFoot,
max_wt_elev_m = as.numeric(max_wt_elev_m))
RealTimeData
folder and annual_appnded_logger_files folder in Boxboglk_all <- rbind(boglk_se_clean, boglk_2019_clean, boglk_2020_clean) %>%
#remove NA values
filter(!is.na(max_wt_elev_m)) %>%
#remove duplicate observations
unique()
#set min and max year
min_year = min(boglk_all$year)
max_year = max(boglk_all$year)
#plot
p_boglk <- ggplot(data = boglk_all) +
geom_line(aes(x = datetime,
y = max_wt_elev_m,
group = 1,
text = paste0("Date: ", date, "\n",
"Max Water Table Elevation (m): ", round(max_wt_elev_m, digits = 3)))) +
theme_classic() +
labs(x = "Time",
y = "Water Table Elevation (m above sea level)",
title = paste0("Bog Lake Fen (BOGLK) Shaft Encoder \n Water Table Elevation (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7))
#interactive plot
ggplotly(p_boglk,
tooltip = "text")
Format the data
#format stripchart data
bog1 <- bogwell_clean_s1 %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart")
#format shaft encoder data
se1 <- s1_all %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#rename columns to join the data
rename(wte = wt_elev_m) %>%
group_by(date) %>%
#calculate max WTE per day
mutate(max_wte = max(wte)) %>%
#remove unnecessary columns
select(-datetime, -wt_elev_ft, -wte) %>%
distinct() %>%
#rename columns to join the data
rename(wte = max_wte)
#join the data
s1_join <- full_join(x = bog1,
y = se1,
by = c("wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
s1_join_sub <- s1_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart))
#set min and max year
min_year = min(s1_join_sub$year)
max_year = max(s1_join_sub$year)
Plot the data
p1 <- ggplot(data = s1_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("S1 1:1 Daily Max WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p1,
tooltip = "text")
summary()s1_diff <- s1_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(s1_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.17488 -0.01869 0.01519 0.09393 0.04801 1.54165
s1_sd <- round(sd(s1_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.35164.
#set min and max year
min_year = min(s1_diff$year)
max_year = max(s1_diff$year)
#plot
ggplot(data = s1_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.002) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("S1 Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
Format the data
#format stripchart data
bog <- bogwell_clean_s2 %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart")
#format shaft encoder data
se <- s2_se_clean %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#remove unnecessary columns
select(-datetime, -max_wt_elev_ft) %>%
#rename columns to join the data
rename(wte = max_wt_elev_m)
#join the data
s2_join <- full_join(x = bog,
y = se,
by = c("wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
s2_join_sub <- s2_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart))
#set min and max year
min_year = min(s2_join_sub$year)
max_year = max(s2_join_sub$year)
Plot the data
p2 <- ggplot(data = s2_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("S2 1:1 WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p2,
tooltip = "text")
summary()s2_diff <- s2_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(s2_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.018694 -0.002111 0.001167 0.001101 0.004255 0.046007
s2_sd <- round(sd(s2_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.00526.
#set min and max year
min_year = min(s2_diff$year)
max_year = max(s2_diff$year)
#plot
ggplot(data = s2_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.002) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("S2 Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
The S3 site was omitted due to unreliable data.
Format the data
#format stripchart data
bog4 <- bogwell_clean_s4 %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart")
#format shaft encoder data
se4 <- s4_2020_clean %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#remove unnecessary columns
select(-datetime, -max_wt_elev_ft) %>%
#rename columns to join the data
rename(wte = max_wt_elev_m)
#join the data
s4_join <- full_join(x = bog4,
y = se4,
by = c("wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
s4_join_sub <- s4_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart))
#set min and max year
min_year = min(s4_join_sub$year)
max_year = max(s4_join_sub$year)
Plot the data
p4 <- ggplot(data = s4_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("S4 1:1 WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p4,
tooltip = "text")
summary()s4_diff <- s4_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(s4_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.0200000 0.0000000 0.0000000 -0.0002602 0.0000000 0.0200000
s4_sd <- round(sd(s4_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.00607.
#set min and max year
min_year = min(s4_diff$year)
max_year = max(s4_diff$year)
#plot
ggplot(data = s4_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.01) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("S4 Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
Format the data
#format stripchart data
bog5 <- bogwell_clean_s5 %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart")
#format shaft encoder data
se5 <- s5_all %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#remove unnecessary columns
select(-datetime, -max_wt_elev_ft) %>%
#rename columns to join the data
rename(wte = max_wt_elev_m)
#join the data
s5_join <- full_join(x = bog5,
y = se5,
by = c("wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
s5_join_sub <- s5_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart))
#set min and max year
min_year = min(s5_join_sub$year)
max_year = max(s5_join_sub$year)
Plot the data
p5 <- ggplot(data = s5_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("S5 1:1 WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p5,
tooltip = "text")
summary()s5_diff <- s5_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(s5_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.0200000 0.0000000 0.0000000 0.0003529 0.0000000 0.0200000
s5_sd <- round(sd(s5_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.00734.
#set min and max year
min_year = min(s5_diff$year)
max_year = max(s5_diff$year)
#plot
ggplot(data = s5_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.01) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("S5 Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
Format the data
#format stripchart data
bog6 <- bogwell_clean_s6 %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart")
#format shaft encoder data
se6 <- s6_all %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#remove unnecessary columns
select(-datetime, -max_wt_elev_ft) %>%
#rename columns to join the data
rename(wte = max_wt_elev_m)
#join the data
s6_join <- full_join(x = bog6,
y = se6,
by = c("wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
s6_join_sub <- s6_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
#unnest the list to turn NULL values to NA values
unnest() %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart)) %>%
#remove duplicate rows
unique()
#set min and max year
min_year = min(s6_join_sub$year)
max_year = max(s6_join_sub$year)
Plot the data
p6 <- ggplot(data = s6_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("S6 1:1 WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p6,
tooltip = "text")
summary()s6_diff <- s6_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(s6_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.080000 -0.010000 0.000000 -0.001182 0.010000 0.030000
s6_sd <- round(sd(s6_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.01252.
#set min and max year
min_year = min(s6_diff$year)
max_year = max(s6_diff$year)
#plot
ggplot(data = s6_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.01) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("S6 Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
Format the data
#format stripchart data
boglk <- bogwell_clean_boglk %>%
#remove unnecessary columns
select(-peatland, -flag) %>%
#create a column to note what type of data this is
mutate(type = "stripchart") %>%
unique()
#format shaft encoder data
se_boglk <- boglk_all %>%
#extract date and year
mutate(date = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y-%m-%d'),
date = as.POSIXct(date, tz = "GMT"),
year = format(as.POSIXct(datetime, format = '%m/%d/%Y %H:%M:%S',
tz = "GMT"),
format = '%Y'),
year = as.numeric(year),
#note data type
type = "shaft_encoder") %>%
#remove unnecessary columns
select(-datetime, -max_wt_elev_ft) %>%
#rename columns to join the data
rename(wte = max_wt_elev_m) %>%
#calculate daily max WTE
group_by(date) %>%
summarise(max_wte = max(wte),
date = date,
year = year,
type = "shaft_encoder"
) %>%
#remove duplicate observations
unique()
#join the data
boglk_join <- full_join(x = boglk,
y = se_boglk,
by = c("wte" = "max_wte", "date", "year", "type"))
#subset the join to the time period where the two data sources overlap
boglk_join_sub <- boglk_join %>%
pivot_wider(names_from = "type",
values_from = "wte") %>%
unnest() %>%
#remove rows with NA values
filter(!is.na(shaft_encoder),
!is.na(stripchart)) %>%
unique()
#set min and max year
min_year = min(boglk_join_sub$year)
max_year = max(boglk_join_sub$year)
Plot the data
p_boglk <- ggplot(data = boglk_join_sub) +
geom_point(aes(x = stripchart,
y = shaft_encoder,
alpha = 0.1,
group = 1,
text = paste0("Date: ", date, "\n",
"Stripchart: ", stripchart, "\n",
"Shaft Encoder: ", round(shaft_encoder, digits = 3)))) +
#plot straight line to compare point spread
geom_abline(intercept = 0,
slope = 1) +
theme_classic() +
labs(x = "Manual Stripchart",
y = "Shaft Encoder",
title = paste0("Bog Lake Fen (BOGLK) 1:1 WTE Comparison (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5,
size = 7),
legend.position = "none")
#interactive plot
ggplotly(p = p_boglk,
tooltip = "text")
summary()boglk_diff <- boglk_join_sub %>%
mutate(diff = shaft_encoder - stripchart)
summary(boglk_diff$diff)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.084811 -0.035978 -0.016262 -0.019127 -0.003954 0.068180
boglk_sd <- round(sd(boglk_diff$diff), digits = 5)
The standard deviation of the differences between stripchart and shaft encoder values is approximately 0.02444.
#set min and max year
min_year = min(boglk_diff$year)
max_year = max(boglk_diff$year)
#plot
ggplot(data = boglk_diff) +
geom_histogram(aes(x = diff),
binwidth = 0.005) +
theme_classic() +
labs(x = "Differences Between Shaft Encoders - Stripcharts (m)",
y = "Frequency",
title = paste0("Bog Lake Fen (BOGLK) \n Daily Max WTE Differences (", min_year, "-", max_year, ")")
) +
theme(plot.title = element_text(hjust = 0.5),
plot.subtitle = element_text(hjust = 0.5))
Note that this file is being knit as index.html
into the climate_bogwell
repository in order to update the GitHub pages website,
where the plots can be viewed online.